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ABSTRACT 


This  report  describes  the  various  techniques  and  algorithms  developed  to 
enhance  the  accuracy  of  the  measurements  of  very  fast  transient  electromagnetic 
fields  generated  during  EMP  testing.  Some  fundamentals  of  digital  signal 
processing  are  presented  and  the  practical  aspects  are  highlighted.  The  response 
of  some  sensors  to  incident  fields  is  studied.  Appropriates  algorithms  for 
reconstructing  the  incident  field  from  the  measurements  are  developed.  The 
effect  of  a  low  frequency  cut-off  of  some  sensors  is  also  studied.  A  method  to 
design  a  digital  filter  of  arbitrary  transfer  function  is  presented.  It  is 
applied  to  compensate  for  losses  occurring  in  long  cables  and  also  to  fine  tune 
other  filters.  The  Hilbert  transform  Is  introduced  to  relate  the  magnitude  and 
phase  component  of  a  transfer  function.  An  implementation  of  the  transform  is 
given,  allowing  the  design  of  digital  filters  when  only  the  magnitude  of  the 
transfer  function  is  known. 


RESUME 


Ce  rapport  decrit  diffdrentes  techniques  developpees  pour  ameliorer  la 
precision  des  mesures  de  champ  dlectromagnetique  extremement  rapide  utilises  lots 
des  test  d' impulsions  41ectromagnetique  (lEM).  Certains  aspects  fondamentaux  du 
filtrage  numerique  sont  revises.  La  reponse  de  capteurs  A  des  champ  electro- 
magnetique  est  etudide.  Differents  algorithmes  de  filtrage  pour  obtenir  le  champ 
incident  A  partir  des  mesures  sont  conqus.  Une  methods  de  conception  de  filtre 
numArique  A  fonction  de  transfert  quelconque  est  presentee.  Cette  methode  est 
appliquee  pour  compenser  1 'attenuation  des  cAbles  et  pour  corriger  les 
imperfections  d'autres  filtres.  La  transformee  de  Hilbert  est  utilise  pour 
etablir  une  relation  entre  1' amplitude  et  la  phase  des  fonctions  de  transfert. 
Cette  relation  permet  de  construire  des  filtres  pour  corriger  des  systAme  dont 
on  ne  peut  mesurer  que  1' amplitude. 
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EXECUTIVE  SUMMARY 


DREO  has  built  a  large  experimental  facility  for  the  measurement  of  the 
effects  of  nuclear  electromagnetic  pulse  (EMP)  on  systems.  The  transient  field 
generated  is  extremely  short  and  difficult  to  measure.  This  report  describes  the 
various  techniques  and  algorithms  developed  to  process  and  enhance  these 
measurements . 

The  response  of  various  sensors  to  incident  fields  is  studied.  It  is  shown 
that  passive  sensors  normally  respond  to  the  time -derivative  of  the  field, 
requiring  some  form  of  integration  to  recover  the  original  field.  Some  well 
known  methods  are  presented.  An  alternative  technique  is  developed,  based  on 
both  hardware  and  software.  It  is  shown  that  the  overall  performances  are 
improved  and  some  experimental  results  are  given.  The  effect  of  a  low  frequency 
cut-off  of  some  sensors  is  also  studied,  and  a  correction  algorithm  is  given. 

A  method  to  design  a  digital  filter  of  arbitrary  transfer  function  is 
presented.  It  is  applied  to  compensate  for  losses  occurring  in  long  cables  and 
also  to  fine  tune  other  filters.  The  Hilbert  transform  is  introduced  to  relate 
the  magnitude  and  phase  component  of  a  transfer  function.  An  implementation  of 
the  transform  is  given,  allowing  the  design  of  digital  filters  when  only  the 
magnitude  is  known. 
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1.0 


Under  a  program  sponsored  by  DND,  DREO  is  expanding  its  Nuclear  Electro¬ 
magnetic  Pulse  (EMP)  research  capability  with  the  addition  of  a  new  large 
experimental  facility.  The  EMP  is  a  high-intensity,  short -duration  electro¬ 
magnetic  field  that  originates  from  a  nuclear  detonation.  Although  the  energy 
content  of  EMP  is  not  very  large  because  of  its  short  duration,  electronic 
components  or  systems  can  be  upset  or  damaged  by  the  very  high  voltages  or 
currents  which  can  be  induced.  To  assure  EMP  hardness,  suitable  tests  and 
analysis  are  needed. 

Even  when  the  EMP  hardening  is  taken  into  consideration  from  the  beginning 
of  the  design  stage  and  supported  by  prediction  and  analysis,  a  validation  of 
such  prediction  by  proper  tests  is  necessary  to  raise  the  level  of  confidence. 
The  EMP  field  or  induced  currents  in  cables  or  antennas  are  simulated  and  the 
system  response  is  compared  to  the  predictions. 

DREO  has  built  an  EMP  test  facility  in  which  a  system  under  test  is 
subjected  to  an  EMP  waveform.  Measurements  are  made  with  a  variety  of  sensors 
and  probes,  and  good  cables  or  fiber-optic  links  are  use  to  carry  these  signals 
into  a  shielded- room  housing  all  the  sensitive  data  acquisition  equipment.  High 
speed  digitizers  must  be  used  to  record  the  extremely  fast  transient  induced  by 
the  EMP  waveform. 

Sensors  have  their  respective  transfer  function,  data  links  introduce 
errors  and  instruments  have  their  limitations,  all  of  which  alter  the  measurement 
being  made.  Post-processing  of  the  data  is  necessary  to  remove  errors  and 
restore  the  original  quantity  that  was  measured  (field  or  current) .  This  report 
describes  the  various  techniques  and  algorithms  developed  to  process  and  enhance 
these  measurements . 

1.1  the  enp  threat 

The  EMP  produced  by  a  nuclear  burst  at  high  altitude  is  a  large  -  ampl i tude , 
very  short  duration  transient  field  covering  a  very  wide  area  beneath  the  burst 
point.  In  general,  the  exact  characteristics  of  the  EMP  field  such  as  peak 
amplitude,  rise  time  and  polarization,  depend  on  many  factors,  such  as  the  weapon 
yield,  the  height  of  burst  and  the  observer's  location.  For  practical  purposes, 
a  standard  waveform  representing  a  worst  case  is  used  [3]  [4]  .  The  standard  EMP 
waveform  has  a  peak  value  of  50  kV  per  meter,  a  rise  time  (10-90%)  of  about 

4.2  nsec,  a  pulse  width  (50-50%)  of  about  185  nsec  and  a  decay  time  (peak-to-10%) 
of  about  600  nsec . 
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Two  analytical  expression  are  conunonly  used  to  approximating  the  EMP 
waveform.  The  double  exponential  form  (1) ,  often  referred  to  as  the  Bell  curve 
[6]  or  the  old  NATO  definition,  is  simpler  and  also  has  a  simple  Laplace 
transform,  which  makes  it  possible  to  solve  simple  problems  analytically: 

E(t)  =  AV  •(e-“*'-e-^'')  (1) 


The  reciprocal  form  (2),  referred  to  as  the  new  NATO  definition,  is  more 
accurate  as  its  derivatives  are  continuous  over  time,  which  gives  a  more 
realistic  leading  edge: 


E(t)  =  AV 


k  at 


l+e 


(2) 


1.2  DREO  LARGE  EKP  SIMDLATOR  FACILITT 

DREO  has  designed  and  built  a  large  bounded  wave  EMP  simulator  of  about  100 
meters  long,  20  meters  wide  and  10  meters  high  with  a  useful  working  area  of  30 
by  20  meters.  Bounded  wave  simulators  produce  electromagnetic  fields  confined 
between  two  parallel  places  or  between  wire  arrays  approximating  two  parallel 
plates.  The  plates  are  tapered  at  both  ends  with  a  generator  attached  at  one  end 
and  a  terminating  load  at  the  other  end.  This  type  of  simulator  is  in  fact  a 
large  transmission  line  with  a  fairly  uniform  field  within  the  parallel  section. 
The  generator  is  capable  of  launching  a  600  kV  impulse  with  a  5  ns  rise  time  in 
the  transmission  line  necessary  to  produce  the  50  kV/m  fields  required  by  the 
MIL-STD-461C  standard. 

The  system  response  is  measured  by  field  sensors  and  current  or  voltage 
probes.  Field  sensors  are  generally  passive  devices,  basically  electric  dipole 
or  magnetic  loop.  Their  response  can  be  very  fast  with  a  high  dynamic  range,  but 
their  sensitivity  is  relatively  low.  Their  output  is  generally  the  time- 
derivative  function  of  the  field  (known  as  D-dot  and  B-dot  sensors)  which  will 
require  further  processing  to  restore  the  waveform  and  for  D-dot  sensors,  the 
output  voltage  can  be  expressed  by: 

V„  =  R.A„..g  =  R.A.,..g.C0Sfl  (3) 

High  speed  digitizers  are  used  to  convert  the  measurements  into  digital 
form.  A  number  of  them  are  used,  at  present  allowing  up  to  a  total  of  nine 
simultaneous  measurements.  The  table  below  summarize  the  key  characteristics  of 
each.  These  instruments  give  the  capability  of  measuring  transients  (single 
shot)  with  a  rise  time  shorter  than  1  or  2  nanoseconds.  Due  to  the  limited 
number  of  samples  obtained  at  the  faster  rates  (>  1  Gs/sec) ,  it  is  not  always 
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possible  to  capture  the  whole  signal,  which  may  last  more  that  1-2  /iS ,  and  still 
capture  accurately  the  fast  rising  edge  (few  nanoseconds).  Two  channels,  at 
different  speed,  may  be  used  to  record  the  whole  signal. 


Model 

Bandwidth 

Sampling  rate 

SCDIOOO 

1  GHz 

200  G- s/sec 

TEK7912HB 

700  MHz 

100  G* s/sec 

DSA602 

400  MHz 

2  G- s/sec 

LeCroy  6880A 

400  MHz 

1.3  G- s/sec 

LeCroy  8828C 

100  MHz 

200  M- s/sec 

The  recorded  data  is  transferred  to  a  computer  for  storage  and  post¬ 
processing.  It  is  very  important  to  save  along  with  the  measurements  all  the 
information  necessary  to  identify  the  setup  used  for  the  test  and  all  the 
instrument  parameters  necessary  to  reconstruct  the  waveform. 

Data  processing  of  the  waveforms  is  an  essential  part  of  the  measurement 
process.  It  is  possible  to  model  arbitrary  transfer  functions  and  to  compensate 
for  errors,  such  as; 

•  recover  from  a  sensor  transfer  function  (by  differentiation, 
integration,  scaling,  etc.) 

•  account  for  sensor  calibration 

•  compensate  for  cable  attenuation  at  higher  frequencies 

•  account  for  any  attenuator  or  filter  inserted 

•  compensate  for  some  of  the  instrument  limitations 

•  merge  two  measurements  of  the  same  signal 

»  perform  functions  on  multiple  measurements  (averaging,  dejitter, 
etc . ) 
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PROCESSIWC 


)  DICITAl.  STCWAL 

Digital  signal  processing  (DSP)  is  the  digital  counterpart  of  analog 
filtering.  Algorithms  to  be  applied  onto  digitized  signals  are  designed  to 
approximated  the  desired  filtering  functions.  Digital  filters  may  be  real  time, 
in  which  case  their  output  is  converted  back  to  analog  signals,  or  they  may 
operate  on  pre-stored  data  sequences  of  finite  length.  The  transfer  function  of 
the  digital  filter  can  be  precisely  tailored  and,  contrary  to  its  analog 
equivalent,  is  not  affected  by  tolerance  or  variation  of  components. 

In  the  discrete- time  linear  systems  theory,  from  which  digital  processing 
is  derived,  signals  are  defined  only  for  the  discrete  values  of  time  t  “  n-T, 
where  T  (also  often  denoted  as  At)  is  the  interval  between  samples.  Discrete¬ 
time  signals  are  represented  as  sequences  of  numbers.  Related  to  the  analog 
world,  a  sequence  may  represent  an  analog  signal  sampled  at  a  regular  interval. 
A  sequence  is  denoted  (x)  and  the  n^**  element  of  it  as  x(n)  or  x„.  Sequences  may 
be  infinite  or  have  a  finite  number  of  samples  N  (xq  to  x^.^) 

The  most  straightforward  approach  to  compute  the  response  of  a  system  is 
to  multiply  the  input  and  the  transfer  function  in  the  frequency  domain: 

e(t)  =»  E(w) 

H(«)  (4) 

U 

y(t)  <=  Y(«) 

and  use  a  fast  Fourier  transform  (FFT)  algorithm  to  generate  the  time  domain 
response.  This  approach  has  several  drawbacks.  e(t)  and  H(w)  must  have  the  same 
size,  which  must  be  a  power  of  2  in  order  to  use  the  standard  FFT  algorithm. 
Although  the  FFT  is  efficient  to  compute  the  spectrum  of  a  signal,  it  is  not  the 
most  efficient  approach  to  filtering,  specially  for  longer  sequences  (some  of  the 
digitizers  may  record  as  many  as  64,000  data  points). 

The  FFT  has  also  some  side-effects ,  which  may  alter  the  solution 
considerably.  Consider  for  instance  a  measurement  truncated  in  order  to  obtain 
details  of  the  rising  edge.  This  truncation  will  introduce  a  severe  aliasing  in 
the  solution.  This  aliasing  is  caused  because  taking  the  inverse  FFT  of  a 
product  corresponds  to  performing  a  circular  convolution  in  time  domain.  This 
problem  can  be  eliminated  if  both  sequences  are  doubled  by  adding  zero's,  which 
would  then  corresponds  to  performing  a  linear  convolution.  This  however  is  at 
the  expense  of  an  increase  of  the  compute  time  and  storage  area  by  a  factor  of  4. 
Nonetheless,  truncating  a  sequence  is  equivalent  to  multiplying  it  with  a 
rectangular  window,  which  corresponds  to  the  convolution  of  H(w)  with  a  sinc(x) 
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function.  This  Is  the  well-known  Gibbs  phenomenon,  which  clauses  H(w)  to  ripple 
and  overshoot  by  as  much  as  9%  around  discontinuities  and  rapid  transitions. 

A  better  approach  to  perform  digital  filtering  is  the  use  of  the 
z-transform,  based  on  the  theory  of  discrete- time  linear  systems. 


2 . 1  THE  Z-TRANSFOBM 

In  continuous -time  systems  theory,  we  use  the  Laplace  and  the  Fourier 
transforms  to  characterize  linear  time -invariant  systems.  Similarly,  it  is 
possible  to  generalized  these  transforms  for  discrete- time  systems,  resulting  in 
the  z- transform. 

The  transform  X(Z)  of  a  sequence  x(n)  is  defined  as; 

a 

X(z)  »  J^x(n)z""  (5) 

where  z  is  a  complex  variable .  The  function  will  converge  only  for  some  values 
of  z.  Detailed  discussion  of  the  z-transform  can  be  found  in  [7],  but  we  can 
state  the  general  results  for  our  cases  of  interest.  The  series  will  converge 
for  values  of  z  on  the  unit  circle  if  the  system  is  stable,  and  vice-versa.  The 
z-transform  will  always  converge  on  the  unit  circle  for  finite-length  sequences 
and  will  converge  for  causal  sequences  if  all  the  poles  of  the  transform  lie 
inside  the  unit  circle. 

The  z-transform  may  be  related  to  the  Fourier  transform  by  evaluating  it 
on  the  unit  circle,  ie.  by  using  the  relation: 

z  *  e^"  (6) 

in  (5)  above,  yielding  to: 


X(z)|x-*>  =  X(e“J)  =  £  x(n)e-j«-"  (7) 

n»-® 

which  is  the  Fourier  transform  of  x(n) .  This  relation  is  important  as  it  gives 
the  exact  frequency  response  of  a  digital  filter  with  a  transfer  function  H(z). 
Similarly,  the  z-transform  may  also  be  related  to  the  Laplace  transform  with  the 
substitution: 


z 


= 


(8) 
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Not  surprisingly,  the  z- transform  shares  the  same  basic  properties  of  the 
Laplace  and  Fourier  transforms.  The  z- transform  is  linear;  that  is  the 
transform  of  axi(n)  +  bx2(n)  is  aXi(z)  +  bX2(z) . 

One  of  the  most  important  and  useful  property  is  the  shift  in  time  (delay) 
of  a  sequence.  If  a  sequence  x(n)  is  delayed  by  no  samples,  then  the  following 
relation  is  obtained: 

x(n-n(,)  z'””X(z)  (9) 

where  z'"“  can  be  seen  as  a  delay  element,  similar  to  e'^°®  with  the  Laplace 
transform  and  e*'^"*^'’  with  the  Fourier  transform. 

As  with  the  Laplace  transform,  it  can  be  shown  that  the  multiplication  of 
two  z- transforms  corresponds  to  the  convolution  of  the  two  corresponding 
sequences ,  ie . : 

x(n)  •  y(n)  X(z)  Y(z)  (10) 

where  the  convolution  is  defined  by: 

x(n)  ♦  y(n)  =  x(k)  y(n-k)  (11) 

k*-« 

Furthermore,  the  input-output  relation  of  a  system  corresponds  to  the 
multiplication  of  the  z- transforms  of  the  input  and  the  unit-sample  response,  in 
a  way  similar  to  the  Fourier  transfoim  as  expressed  by  (4) .  This  property 
suggest  that  one  could  use  the  same  method  for  solving  problems.  It  would 
require  performing  the  z- transform  of  the  input  with  (5)  and  an  inverse  transform 
of  the  output  (not  shown,  but  it  is  a  rather  complicated  function,  difficult  to 
calculate)  .  Instead,  if  we  define  H(z)  as  a  ratio  of  two  polynomials  in  z*^  such 
as: 
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which  may  be  normalized  for  ao-1  and  rewritten  as^: 


"  boX„  +  biX„.i  +  bjXn.j  +  •  •  ■  bMXn-M  (14) 

■  ®iyn-l  ~  “  •  •  •  %yn-N 

which  states  that  the  output  at  a  given  time  n  is  a  function  of  the  present  and 
past  values  of  the  input,  and  also  a  function  of  the  past  values  of  the  output. 
The  problem  of  designing  a  filter  of  given  characteristics  becomes  a  question  of 
finding  the  proper  coefficients. 

All  the  algorithms  developed  in  this  report  are  based  on  the  relations 
between  the  z- transform  and  the  Fourier  and  Laplace  transforms  and  the  use  of 
several  techniques  to  transpose  a  transfer  function  from  the  frequency  domain  or 
the  s-domaln  into  the  z-domain  to  obtain  H(z).  Relation  (14)  provides  an 
efficient  computer  implementation  of  the  resulting  filter. 

Several  techniques  have  been  developed  to  design  digital  filters,  yielding 
to  two  general  classes  of  filters.  For  filters  where  A(z)  is  1  (ag  is  1,  all 
other  terms  are  0) ,  the  output  depends  on  the  input  only  and  consequently  the 
impulse  response  of  such  system  will  be  of  finite  duration.  Such  filterr  are 
known  as  Finite  Impulse  Response  (FIR)  filters.  For  cases  where  the  output  may 
also  depend  on  the  past  values  of  the  output,  the  impulse  response  is  infinite 
because  of  the  recursive  nature  of  these  filters.  'These  filters  are  known  as 
Infinite  Impulse  Response  (HR)  filters.  FIR  and  HR  filters  each  have  their 
advantages  and  disadvantages  and  the  choice  between  one  or  another  and  the 
technique  used  to  find  the  coefficients  is  often  dictated  by  the  particular 
problem: 

•  HR  filters  are  the  most  efficient,  requiring  5  to  10  times  fewer 
coefficients  than  the  FIR  equivalent. 

•  HR  filters  are  primarily  limited  to  more  classical  filters,  such  as 
low-pass,  high-pass,  pass-band,  band-stop,  etc.,  for  which  the 
expression  in  Laplace  domain  is  known.  It  also  allows  the  use  of  the 
classical  solutions  of  filter  design,  such  as  Butterworth, 
Chebyshev,  elliptic  filters. 

•  An  expression  for  calculating  the  coefficients  of  HR  filters  for 
any  AT  can  usually  be  obtained.  By  comparison,  the  coefficients  of 
FIR  filters  are  often  obtained  by  iterative  or  cumbersome  procedures 
and  are  calculated  for  a  single  AT  only,  therefore  requiring  the 


^  A  detailed  development  of  this  relation  can  be  found  in  [7]. 
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storage  in  advance  of  tables  containing  the  coefficients  of  various 
AT. 

FIR  filters  can  be  designed  to  have  a  precisely  linear  phase,  which 
is  not  possible  with  HR  filters. 

FIR  filters  can  be  designed  to  match  any  arbitrary  response. 


2.2  DESIGN  OF  HR  FILTERS 

Designing  a  digital  filter  (ie.  finding  the  proper  coefticients)  is 
basically  a  problem  of  approximating  the  characteristics  of  the  desired  filter 
by  following  some  criteria.  HR  filters  are  best  suited  for  cases  where  an 
equivalent  analog  filter  exists.  Several  techniques  (such  as  the  impulse 
invariance,  mapping  of  differential  equations,  the  bilinear  transformation,  and 
the  matched  z  transformation)  have  been  developed  to  transpose  a  filter  defined 
in  s-domain  into  the  z-domaln.  Of  particular  interest  to  our  application  is  the 
bilinear  transform. 


This  procedure  is  developed  in  [7]  by  solving  a  simple  problem  with  the 
integral  equation  and  approximating  the  integral  with  the  trapezoidal  rule, 
resulting  in  the  substitution^: 


s 


2  1-z-i 


(15) 


The  main  advantages  of  the  bilinear  substitutions  are; 

•  A  order  transfer  function  H(s)  gives  a  N***  order  digital  filter, 
thus  simple  filters  will  result  into  simple  equations. 

•  A  stable  analog  filter  will  always  give  a  stable  digital  filter. 

•  The  transformation  does  not  Introduce  aliasing.  Therefore,  it  is 
possible  to  design  filters  operating  near  their  Nyquist  limit. 

There  is  however  one  disadvantage,  that  is  the  introduction  of  a  distortion 
of  the  frequency  axis.  The  relation  between  the  analog  frequency  0  and  the 
digital  frequency  w  is  no  longer  linear,  but  corresponds  to: 


^  This  expression  is  also  an  approximation  of  (8),  ie.  the  truncation  to  one 
term  of  the  Taylor  series  expansion  of  s  *  Inz  /  T. 
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OT 


(16) 


which  may  be  considered  linear  for  low  frequencies  and  the  distortion  will  grow 
as  the  frequency  approaches  the  Nyquist  frequency.  For  instance,  for  1/10, 
1/4  and  1/2  of  the  Nyquist  frequency,  the  distortion  is  1%,  5%  and  15% 
respectively.  In  many  cases,  this  distortion  may  be  accounted  for.  Consider  for 
Instance  the  design  of  a  low-pass  or  band-pass  filter  with  a  critical  cut-off 
frequency  It  is  possible  to  still  obtain  a  very  acceptable  filter  with  the 
bilinear  transformation  by  solving  the  problem  for  the  cut-off  frequency 
(2/T)tan(WcT/2)  instead  of  w,,,  ie.  by  pre-distorting  the  frequency  axis. 


2.3  DESIGN  OF  FIR  FILTERS 

FIR  filters  are  attractive  particularly  for  their  capability  to  implement 
arbitrary  frequency  responses.  As  seen,  the  transfer  function  of  a  FIR  filter 
is  given  as; 


M-l 

H(z)  ■  J^h(n)2'”  (17) 

n-O 

where  h(n)  is  the  Impulse  response  of  the  filter.  A  FIR  filter  is  designed 
simply  by  finding  its  impulse  response^.  However,  in  most  cases,  the  impulse 
response  is  infinite,  ie.  it  extends  from  -w  to  +«.  The  most  straightforward 
approach  consists  of  truncating  the  impulse  response  to  N  terms.  This 
truncation,  however,  leads  to  the  well-known  Gibbs  phenomenon.  Simply, 
truncation  in  time  domain  is  similar  to  a  multiplication  with  a  rectangular 
window  and  results  in  frequency  domain  in  the  convolution  of  H(z)  with  the 
transform  of  the  window,  which  is  sin(wN/2)/sin(w/2) .  The  result  is  a  smear  of 
the  frequency  response,  especially  near  the  sharp  transition,  resulting  in  a 
fixed  percentage  of  overshoot  and  ripple.  For  example,  the  overshoot  of  an  ideal 
low-pass  filter  is  about  9%  at  the  transition,  and  more  important,  it  does  not 
decrease  as  the  size  of  h(t)  is  increased.  The  problem  may  be  controlled  with 
the  use  of  other  windows,  such  as  the  Hamming,  Hanning,  Blackman,  Bartlett  and 
Kaiser  windows.  These  windows  increase  the  stop-band  attenuation  and  decrease 
the  main  lobe  width,  resulting  in  smaller  overshoot  and  a  more  abrupt  stop  band. 
This,  however,  is  at  the  cost  of  a  larger  N.  In  our  applications,  where  we  are 
concerned  a'jout  compensating  for  cable  losses,  etc.,  the  frequency  response  is 
fairly  smooth  and  thus,  the  Gibbs  effect  is  very  small  and  may  be  neglected. 


^  Another  technique,  known  as  the  frequency- sampling,  allows  to  design  FIR 
filters  by  specifying  N  samples  of  their  frequency  response.  This  design 
is  particularly  attractive  for  narrow  pass -band  filters,  but  in  general 
the  windowing  techniques  are  more  practical  and  more  flexible . 
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Ideally,  the  analytical  form  of  h(n)  should  be  used  to  calculate  the 
coefficients.  In  many  cases,  it  is  difficult  to  evaluate  and  a  M-point  inverse 
DFT  or  FFT  is  used  with  a  N-point  window.  This  gives  accurate  results, 
especially  if  the  frequency  response  is  smooth  and  if  M  »  N.  However,  great 
care  should  be  exercise  to  avoid  the  side-effects  resulting  from  the  application 
of  the  DFT  algorithm.  In  any  case,  relation  (6)  should  be  used  to  check  the 
response  of  the  filter,  specially  for  frequencies  not  a  multiple  of  Af  as  the 
filter  may  give  an  exact  solution  for  these  frequencies  but  exhibit  some  ripple 
for  others.  A  technique  for  designing  FIR  filters  will  be  further  developed  in 
Chapter  4. 

Some  problems  are  more  difficult  (or  impractical)  to  implement  as  FIR 
filters  because  their  impulse  response  is  of  long  duration  and  would  require  a 
large  number  of  points.  Consider  for  instance  a  1*^  order  low-pass  filter  with 
a  cut-off  frequency  Wg-O.Ol/T.  Its  impulse  response  is  exp(-a),.nT)  and  between 
250  and  500  samples  would  be  required  to  create  a  FIR  filter.  Filters  which 
exhibit  resonances  (narrow  band  filters)  also  have  a  long  duration  impulse 
response.  On  the  other  hand,  filters  with  a  near  unity  frequency  response,  such 
as  one  designed  to  compensate  for  the  slight  loss  occurring  in  cables  at  higher 
frequencies,  have  a  very  short  impulse  response^,  resulting  in  typically  between 
5  to  20  coefficients. 

This  last  property  can  be  used  to  combine  both  FIR  and  HR  techniques  to 
design  more  efficient  filters.  Consider  for  instance  the  problem  shown  on 
Figure  2-1  where  the  desired  function  H(w)  cannot  be  expressed  easily  in  the 
s- domain.  A  FIR  design  would  give  a  very  good  approximation  of  the  function,  but 
with  a  large  number  of  coefficients.  Instead,  H(w)  can  be  represented  as  the 
product  two  functions:  G(w)  being  an  ideal  1®^  order  low-pass  filter,  and 
E(a))  -  H(w)/G(w)  being  the  error  between  the  two.  G(w)  is  efficiently 
implemented  as  a  1®^  order  HR  filter  and  E(w)  as  a  FIR  filter  with  few 
coefficients.  This  results  obviously  in  an  HR  filter. 


^  The  impulse  response  of  a  pass-all  filter,  H(w)  -  1,  is  the  ideal  impulse 
5(n).  Consequently,  the  impulse  response  will  tend  toward  the  ideal 
Impulse  as  the  frequency  response  approaches  unity. 
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3.0  CORBECTIWG  FOR  ELECTROMAGNETIC  FIELD  MEASTOtEMEHTS 


Field  sensors  are  used  to  measure  the  electric  field  (E,  in  Volt/meter)  or 
flux  density  (D-cqE,  in  Coulomb/meter ,  where  Cq-S . 85 •  10'^^)  and  the  magnetic  field 
(H,  in  Ampere/meter)  or  flux  density  in  Tesla,  where  /io“4w •  10'^) .  Field 

sensors  may  be  active  or  passive  (ie.  no  external  power),  but  generally  passive 
devices  are  used.  The  most  basic  type  of  field  sensor  is  the  small^  electric 
dipole  or  the  small  magnetic  dipole  (loop). 

Figure  3-1  shows  a  simple  circuit  model  for  a  small  electric  and  a  small 
magnetic  dipoles.  The  transfer  function  of  the  electric  dipole  is  given  as: 

^out  "  tiaq*  ° 

which  both  output  a  voltage  proportional  to  the  normal  component  of  the  field  (d 
being  the  angle  between  the  field  and  vector  normal  to  the  sensor  ground  plane) 
above  a  certain  cut-off  frequency  (related  to  the  time  constants  RC  and  L/R) ,  and 
a  voltage  proportional  to  the  time -derivative  of  the  field  below  it.  For 
practical  geometries  and  for  sensors  normally  matched  to  the  instrumentation 
(50  0) ,  these  sensors  operate  well  below  the  cut-off  frequency  and  thus  will 
always  respond  to  the  derivative  of  the  field,  approximated  by: 

Vout  “  A,^*^cos«  (19) 


Similarly,  the  output  of  the  magnetic  dipole,  is  given  by: 

* 

Because  these  sensors  output  the  time -derivative  of  the  excitation,  they 
are  known  as  D-dot  and  B-dot  sensors  (or  D  and  B  sensors).  Since  the  response 
of  these  sensors  depend  on  the  geometry  and  the  external  load  only,  accurate 
calibration  can  be  obtained  from  the  sensor  dimensions  solely. 

It  is  obviously  necessary  to  perform  an  integration  on  D-dot  and  B-dot 
sensors  output  to  restore  the  original  field  waveform.  Section  3.1  will  present 
several  methods  to  achieve  Integration.  Section  3.2  will  discuss  the  case  of  a 


=  A, 


eq - T — 


>B  cos^ 


(20) 


^  Small  is  define  as  electrically  small,  relative  to  the  wavelength  of 
highest  frequency.  For  standard  EMP  waveform,  it  corresponds  to  sensors 
smaller  than  10-15  cm. 
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Electric  dipole  model  Magnetic  dipole  model 


Figure  3-1.  Simple  models  of  electromagnetic  field  sensors. 
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small  short-circuited  loop  used  to  measure  the  H  field  directly. 


3.1  niTEGRATIIIG  D-DOT  AND  B-DOT  SENSORS  ODTFOT 

Two  traditional  and  straightforward  methods  have  been  used  to  integrate  the 
measurements  of  D-dot  sensors^:  using  software  integration  and  using  hardware 
integration.  A  better  method  combining  the  advantages  of  both  will  be  presented 
in  this  section. 

This  simplest  of  all  methods  is  to  digitize  the  raw  data  from  the  sensor 
and  perform  the  integration  by  software  with  the  simple  algorithm: 


Yn  *  Yn-i  + 


(21) 


This  technique  works  well,  but  it  has  major  disadvantages.  Because  the 
D-dot  sensors  respond  to  the  time - derivative ,  they  produce  an  output  which  is 
even  faster  that  the  excitation.  Figure  3-2  (top)  shows  in  (a)  a  typical  D-dot 
measurement  with  a  rise  time  of  approximately  1.2  ns,  which  is  about  5  times 
faster  than  the  actual  field  being  measured  as  shown  on  bottom^.  It  is 
therefore  necessary  to  oversample  the  signal  to  capture  the  faster  rising  edge. 
An  other  major  problem  with  this  technique  is  the  reduction  of  the  dynamic  range 
resulting  in  a  poorer  resolution  of  the  part  of  the  signal  corresponding  to  the 
falling  edge  of  the  pulse  (this  problem  will  be  discussed  later) .  This  algorithm 
is  also  highly  susceptible  to  the  trigger  jitter  unless  a  sufficient  number  of 
points  is  obtained  on  the  rising  edge.  Of  course,  this  technique  also  require 
higher  bandwidth  instrumentation  and  data  links,  which  becomes  very  expensive  at 
frequencies  above  1  GHz. 

The  other  common  technique  to  integrate  the  measurement,  is  the  use  of  a 
passive  integrator.  They  consist  of  simple  resistor-capacitor  passive  network 
resulting  in  a  l**^-order  low-pass  filter  with  the  transfer  function: 

H(s)  -  ^ —  (22) 

1  +  r  s 

which  approximate  an  integrator  for  frequencies  large  compared  to  1/r .  For 
commonly  available  devices  (r-  1  to  10  /isec) ,  the  rising  edge  will  be  properly 
restored  and  consequently,  there  is  no  need  to  oversample  as  the  signal  being 


This  section  also  applies  to  B-dot  magnetic  sensors. 

^  Actually,  the  reference  curve  shown  is  the  voltage  measured  across  the  two 
plates  of  the  simulator  near  the  generator.  This  voltage  is  faster  than 
the  actual  field  found  at  the  sensor  location,  which  is  about  6.5  ns. 
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measured  is  not  much  faster  than  the  original  field.  On  the  other  hand,  the 
falling  edge  will  sag  due  to  the  cut-off  frequency  in  the  range  of  30-300  kHz, 
and  will  need  to  be  restored  with  further  processing.  A  disadvantage  of  these 
devices  is  that  they  are  designed  to  match  high- impedance  (1  MO)  instrumentation, 
which  is  not  usual  for  high-bandwidth  digitizers^.  But  the  biggest  disadvantage 
is  their  high  insertion  loss.  The  output  of  a  D-dot  sensor  to  a  40  kV/m  field 
is  about  25  Volt  peak,  and  the  output  of  a  1  /xsec  integrator  is  75  mV  peak.  It 
then  becomes  impossible  to  measure  lower  fields,  such  as  in  boxes  or  enclosures. 

An  alternate  technique  developed  to  overcome  these  disadvantages  make  use 
a  specially  designed  device  similar  to  the  passive  integrator,  but  with  a  much 
smaller  time  constant  of  the  order  of  2  to  10  ns,  that  has  been  named  a  partial 
integrator.  The  result  is  that  only  the  higher  frequencies  which  contribute  for 
the  accelerated  rise  time  are  integrated.  The  output  still  needs  to  be 
processed,  but  it  is  only  10  to  50%  faster  that  the  original  field,  making  it 
easer  to  digitise.  The  effect  of  such  a  device  can  be  observed  on  Figure  3-2 
(top),  where  the  curve  (b)  is  the  output  of  a  5  ns  partial  integrator.  The 
sensor  output  is  slowed  to  give  a  rise  time  of  about  4.1  ns,  thus  allowing  the 
use  of  a  slower  sweep  to  capture  a  longer  time  window. 

Another  advantage  of  the  passive  integrator  is  the  increase  in  dynamic 
range.  The  time -derivative  of  the  field  result  in  a  signal  of  large  amplitude 
during  the  short  rising  portion  of  the  field,  followed  by  a  much  smaller  signal 
(in  this  example,  with  an  amplitude  of  1/25^**  of  the  peak)  of  long  duration 
during  the  decay  portion.  The  digitizers  have  a  resolution  of  8  to  9  bits  for 
a  full  scale  signal  necessary  to  capture  the  peak  amplitude,  leaving  a  resolution 
of  only  about  4  bits  for  the  trailing  edge.  This  poor  resolution  will  introduce 
noticeable  errors  during  the  integration.  In  this  example,  the  5  ns  partial 
integrator  reduces  the  peak  amplitude  by  only  40%,  but  a  30  ns  partial  integrator 
would  reduce  the  amplitude  by  a  factor  of  5,  resulting  in  an  increase  in  dynamic 
range  of  about  2  or  3  bits. 

To  build  even  a  simple  filter  for  frequencies  above  100  MHz  is  difficult. 
Discrete  components  cannot  be  used  because  of  stray  capacitance  of  inductors  and 
stray  inductance  of  leads  and  capacitors.  Even  very  short  leads  introduce  stub 
effect.  This  problem  has  been  addressed  by  designing  a  transmission  line  section 
about  5  inches  long.  The  diameter  of  the  central  conductor  and  the  outer  shell 
are  about  0.312  and  0.718  inches  respectively  and  are  chosen  to  maintain  the  50  n 
characteristic  impedance.  The  middle  portion  of  the  central  conductor  is  made 
much  thicker,  almost  touching  the  outer  shell,  thus  inserting  a  relatively  large 
capacitance.  This  conductor  is  made  of  aluminum  which  is  anodized  to  form  a  very 


^  The  LeCroy  digitizers  have  a  front-end  amplifier,  which  may  be  used  to 
provide  a  high  Impedance  input,  but  at  the  cost  of  limiting  the  bandwidth 
to  150  MHz. 
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good  Insulating  layer  at  the  surface  to  prevent  any  conduction  to  the  outer 
shell.  This  capacitance  (C^)  appears  in  parallel  with  the  load  and  values  in  the 
range  of  100  to  400  pF  can  be  obtained.  This  capacitance  is  much  larger  than  the 
internal  capacitance  C  of  the  sensor,  changing  its  transfer  function  to: 

Vout  =  A^R  cosfl  .  for  C,  »  C  (23) 

which  is  equivalent  to  the  sensor  output  passed  through  a  l®*^-order  filter  with 
a  cut-off  frequency  f^  of  l/2wRC,.  It  is  difficult  to  build  the  device  to  match 
a  given  f^,  but  as  long  as  the  transfer  function  can  be  measure  accurately, 
proper  waveform  recovery  can  be  achieved.  Figure  3-3  shows  the  measured  transfer 
function  of  a  partial  integrator  compared  with  an  ideal  1®*^- order  low- pass 
filter,  with  at  cut-off  frequency  of  30  MHz.  The  agreement  is  very  good:  less 
than  1  dB  deviation  for  frequencies  up  to  750  MHz,  making  it  suitable  for  an  HR 
implementation^ . 

Recovery  of  the  original  field  is  done  simply  by  multiplying  by 
(1  +  RCjS  )  to  undo  the  effect  of  the  partial  integrator  (Cjj  «  100  pF)  ,  and  then 
multiply  by  1/s  to  integrate.  By  using  the  bilinear  transformation,  we  obtain 
the  transfer  function: 


H(z) 


(|+RC,)  +  (|-RC+x)z*^ 

i  -  z 


(24) 


which  translates  to: 

yn  -  yn-i  *  -  (-^-Rc,)x„-i  (25) 

On  Figure  3-2  (bottom),  the  result  of  this  algorithm  is  shown  on  curve  (b). 
The  agreement  with  a  purely  software  integration  on  curve  (a)  is  very  good.  The 
algorithm  is  very  stable,  as  it  can  be  seen  on  curve  (c) ,  where  the  sweep  rate 
was  selected  to  be  very  near  the  Nyquist  limit  (sampling  at  about  every  2  ns) . 
The  key  characteristics  of  the  field  (rise  time,  peak  amplitude)  and  the  overall 
shape  of  the  curve  is  still  very  close  to  the  other  measurements  (done  at  a 
0.2  ns  sampling  rate),  even  when  severe  blooming  occur  onto  the  internal 
graticule  of  the  digitizer  such  as  in  this  measurement. 


^  The  small  residual  error  can  be  easily  corrected  with  a  small  FIR  filter 
if  necessary. 
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3.2  amPEHSATING  FOR  IjOV  FSEQUEHCT  CUT-OFF 

It  was  shown  that  passive  sensors  connected  to  50  fl  lines  gives  the  time- 
derivative  of  the  fields.  However,  if  we  consider  the  case  of  a  small  short- 
circuited  loop,  it  is  possible  to  obtain  magnetic  fields  directly  by  measuring 
the  current.  Equation  (20)  then  becomes  ==  Aaq*B/L.  This  response  is  no 
longer  dependent  of  the  geometry  alone  and  the  sensor  must  be  properly 
calibrated.  The  main  parameter  of  the  loop,  beside  its  equivalent  area,  is  its 
inductance,  which  is  given  by; 

L  -  Morln(8-|-2)  (26) 

which  for  small  loops  is  of  the  order  of  0.1  /xH.  The  DC  resistance  of  the  loop 
is  negligible  compared  to  the  impedance  of  the  probe.  A  current  probe  is 
basically  a  transformer  whose  primary  is  the  wire  carrying  the  current  to 
measure.  Typical  current  probes  have  an  insertion  impedance  of  20  to  50  mfl. 

The  values  of  R  and  L  will  effectively  impose  a  low  frequency  cut-off  at 
fe  -  R/2wL,  in  the  range  of  30  to  300  kHz.  The  effect  of  this  cut-off  can  be 
seen  on  Figure  3-4.  The  excitation  (reference  curve)  is  an  ideal  EMP  waveform, 
with  a  rise-time  of  6.5  ns,  a  pulse  width  of  266  ns  and  duration  of  about  1  ixs , 
matching  the  EMP  simulator  characteristics.  Not  surprisingly,  the  sensor  will 
give  an  accurate  measurement  of  the  rising  edge,  but  the  falling  edge  will  be 
affec  .>d  by  the  cut-off.  The  waveform  crosses  the  zero-axis  at  a  earlier  time 
and  converges  back  to  zero  from  the  negative  side. 

Manufacturers  of  field  sensors  suggest  a  method  to  correct  this  distortion 
[11],  but  it  is  derived  from  integrals  in  their  continuous  form,  which  is  not 
very  practical.  Instead,  the  correction  may  be  achieved  by  using  the  HR  filter 
equivalent  of  (s+2»rfg)/s,  yielding  to: 

Yn  *  Yn-i  +  (2xfc-2*l)  x„  -  (2jrf^.J-l)  x„.i  (27) 

Numerical  error  may  be  introduced  for  lower  value  of  f^  or  T  because  the 
coefficient  of  x„  and  x„-i  then  become  very  close  to  unity.  This  problem  can  be 
overcome  by  using  double  precision  arithmetic. 

As  noted,  it  is  necessary  to  calibrate  the  sensor,  in  particular  to  find 
the  value  of  f^.  It  is  not  always  possible  to  measure  the  sensor  response  for 
very  low  frequencies^.  It  is  however  possible  to  obtain  this  value  from  a  time 


^  For  instance,  the  HP8753  network  analyzer  used  for  sensor  calibration 
covers  the  band  from  300  kHz  to  6  GHz. 
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Sensor  CXitput 


Figure  3-4. 


Effect  of  the  low  frequency  cut-off  of  a  sensor. 
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domain  measurement.  As  seen  on  Figure  3-4,  the  sensor  output  crosses  zero  at  a 
time  to  and  this  value  is  related  to  f^.  For  to  much  larger  than  the  rise  time, 
the  reference  curve  may  be  approximated  by  e'“*’  (from  equation  (1),  where  p»a) 
and  the  value  of  a  can  be  obtained  from  the  curve.  Solving  the  Laplace  equation 
yield  to: 


a  e  =  Wj  e 

which  can  be  solved  numerically  by  an  iterative  procedure . 


(28) 


22 


4.0  CORRECTHIG  FOR  MEASURKMEHT  SYSTEM  ERRORS 


Beside  specific  transfer  functions  and  errors  introduced  by  the  sensors, 
the  other  components  used  in  making  the  measurement  system  (cables,  fiber  optic 
links,  digitizers,  etc.)  also  introduce  errors.  Many  of  these  errors  can  be 
precisely  measured  and  an  appropriate  filter  to  compensate  their  effect  can  be 
des.j.^-'ed. 

The  first  section  addresses  the  problem  of  losses  in  long  cables.  A 
technique  is  developed  to  design  FIR  filters  of  arbitrary  transfer  functions. 
Although  the  problem  with  cables  is  not  severe,  the  same  technique  is  used 
whenever  a  precise  filter  is  required,  such  as  in  the  case  to  fine  tune  an  HR 
filter. 

The  second  section  treats  the  problem  of  designing  a  digital  filter  when 
only  the  desired  magnitude  is  known.  The  Hilbert  transform  is  introduced  to 
relate  the  magnitude  and  the  phase  response  of  a  system. 


4.1  COSBECTING  FOR  CABLE  LOSSES 

The  DREO  large  EMP  simulator  has  a  test  area  of  30  by  20  meters  located 
about  15  m  from  the  shielded  room  housing  the  instrumentation.  Therefore,  analog 
data  links  (coaxial  cables  or  fiber  optic)  of  up  to  40-50  meters  may  be  needed. 
Even  good  quality  cables  fitted  with  good  quality  connectors  attenuate  the  higher 
frequencies  by  as  much  as  3  to  6  dB.  This  section  describes  the  design  of  a 
corrective  filter  to  compensate  for  this  effect. 

Figure  4-1  shows  the  frequency  response,  H(w) ,  of  a  30  m  cable.  The 
apparently  erratic  behaviour  of  the  phase  is  in  fact  a  nearly  linear  phase,  but 
with  a  rapid  slope  due  to  the  long  delay  introduced  by  30  m  of  cable^.  This 
delay  does  not  affect  the  solution  and  it  can  be  removed  easily  (see  below).  The 
inverse  of  this  response,  H5(w)  -  1/H(w) ,  is  used  to  compute  the  recovery  filter. 
In  this  case,  a  FIR  design  was  chosen  to  build  the  filter  as  it  allows  to  model 
arbitrary  transfer  functions.  The  procedure  is  based  on  truncating  the  impulse 
response,  h(n) ,  obtained  from  the  inverse  FFT  if  the  frequency  response.  The 
truncation  is  performed  by  keeping  all  the  points  around  the  peak  until  h(n) 
permanently  stay  below  a  predetermined  threshold,  usually  0.1%  to  1%  of  the  peak. 
This  procedure  is  valid  if  the  sequence  obtained  is  much  shorter  than  the  size 
of  the  FFT,  which  is  the  case  if  the  filter  results  in  fewer  than  50  coefficients 
and  a  FFT  size  of  512  to  2048  is  used.  Because  the  frequency  response  is  very 


^  In  fact,  when  measuring  the  response  of  long  cables  with  the  HP8753 
network  analyzer,  a  very  slow  sweep  should  be  chosen  to  remove  the 
transient  effect  of  the  long  delay. 
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smooth,  the  Gibbs  phenomenon  is  negligible  and  the  use  of  a  window  function  is 
not  necessary. 

It  is  known  that  the  straight  application  of  the  DFT  way  gives  erroneous 
results.  This  subject  has  been  treated  extensively  in  the  literature  (such  as 
[9],  [12],  [13],  [14]  and  [17])  and  will  not  be  reviewed  here. 

The  Fourier  transform  of  a  real  signal  yield  to  a  symmetry  between  the 
positive  and  negative  frequencies;  G(f)  -  G*{-f),  Consequently,  it  is  necessary 
to  supply  both  positive  and  negative  frequencies  to  recover  a  real  signal  from 
its  spectrum.  The  DFT  of  a  N-point  sequence,  from  t-0  to  (N-l)At,  results  in  a 
sequence  of  Va+1  frequencies,  from  f-0  to  followed  by  its  conjugate 
image,  Vz*!  frequencies,  from  f-- (“/2*l-)^f  to  -Af  It  is  therefore  necessary 
to  mirror  the  frequency  response  with  its  conjugate  before  taking  its  IDFT  to 
obtain  a  real  result.  In  addition,  it  is  necessary  that  the  response  for  f— 0  and 
f“(V2)^f  be  real  (these  are  the  only  two  frequencies  which  do  not  have  their 
conjugate  mirrored) .  The  value  at  DC  is  always  real  and  can  be  obtained  from  DC 
measurement  or  obtained  by  extrapolating  the  first  few  frequencies  to  zero.  The 
value  at  f-Vz^f  however  is  generally  complex.  One  possible  work-around  is  to 
use  a  window  function  such  as  a  cosine-taper  to  force  the  frequency  response  to 
zero  or  unity  at  that  frequency.  An  other  solution  is  to  multiply  the  frequency 

response  with  the  ideal  delay  line  e*'^“‘'’.  This  function  has  a  unit  magnitude 
and  a  linear  phase,  corresponding  to  a  shift  of  the  time  domain  respond  by  t-ttj. 
Choosing  t<j  to  obtain  a  response  of  -1  at  f-Vz^f"  ensures  that  the  IDFT  yields 
to  a  real  result,  and  furthermore  ensures  that  the  phase  will  be  continuous,  ie. 
will  not  jump  from  a  negative  to  a  positive  phase  at  that  frequency. 

Of  course,  the  Impulse  response  h(n)  must  be  computed  with  the  same 
sampling  rate  as  the  sequences  to  be  filtered,  which  is  determined  by  the 
specific  digitizer  used.  For  the  HB7912  or  SCDIOOO  for  instance,  the  sampling 
rate  is  related  to  the  time  base  as:  AT  —  lOtj^/N^,  where  is  the  sweep  rate 
in  second  per  division  and  Nt  is  the  number  of  points  (512  or  1024)  .  The 
corresponding  Af  required  is:  Af  -  l/(NfAT)  -  N^/Nj  /10t,„,  where  Nf  is  the  size 
of  the  FFT.  An  appropriate  value  for  Af  is  5  MHz,  for  sweep  rates  of  the  order 
of  10  nsec/div  It  is  not  necessary  to  measure  the  cable  characteristics  for  all 
possible  values  of  Af;  instead,  decimation  or  interpolation  of  only  few 
measurements  may  be  done . 

Figure  4-2  shows  the  function  to  compensate  for  the  cable  losses.  A 
typical  measurement  with  the  HP8753  network  analyzer  is  made  at  401  frequencies, 
and  with  a  Af  of  5  MHz,  it  covers  a  band  up  to  2  GHz.  In  this  example,  the 
coefficients  are  computed  for  a  10  ns/div  sweep  (AT  -  0.195  ns)  and  for  a  FFT 


1 


The  DFT  of  a  sequence 
Therefore ,  the  negative 
(Vz+DAf  to  (N-l)Af. 


is  periodic  itself, 
frequencies  -(Vz'l)^^ 


with  a  period  of  1/AT. 
to  -Af  are  reflected  to 
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Frequency  response  of  a  30  meter  cable . 
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size  of  1024.  Therefore,  it  is  necessary  to  specify  H(,(f)  up  to  2.5  GHz  and  some 
extrapolation  is  required.  One  solution  is  to  use  a  window  function  and  force 
to  converge  toward  zero,  let  say  in  the  band  of  1.5  to  2  GHz.  This  is 
adequate  if  we  consider  that  most  of  the  signal  is  concentrated  below  1  GHz  and 
anything  above  is  considered  as  noise.  Another  solution  is  to  extrapolate  H^Cf) 
either  numerically  (linear  or  quadratic)  or  by  modelling  it  as  a  1®*’  or  2"'*- order 
system.  However,  this  method  may  produce  unstable  filters  as  the  noise  present 
in  the  higher  frequencies  will  be  amplified  (correcting  for  losses  of  3  to  6  dB 
generally  produce  good  results,  but  great  caution  should  be  exercised  for 
correction  of  more  than  10  dB) .  The  method  used  in  this  example ,  is  to  extend 
H5(f)  based  on  the  value  of  the  last  frequencies  present  in  the  measurement.  Its 
overall  performances  are  good  as  the  amount  of  boosting  is  limited  and  it  tends 
to  produce  filters  with  fewer  coefficients. 

An  impulse  response  obtained  form  the  IDFT  is  also  periodic,  ie.  a  circular 
sequence  is  obtained.  Due  to  inherent  or  added  delays,  the  origin  of  the 
response  is  not  necessarily  the  beginning  or  the  middle  of  the  sequence .  The 
easiest  procedure  is  to  perform  a  circular  shift  to  bring  the  peak  value  near  the 
middle,  then  perform  the  truncation  and  possibly  multiply  with  a  window  function. 
Figure  4-3  shows  the  function  h(t)  obtained  from  Hc(f)  on  Figure  4-2  (only  about 
80  coefficients  are  shown).  By  using  a  threshold  of  1%  of  the  peak  value,  this 
sequence  can  be  truncated  to  9  coefficients.  The  exact  frequency  response  of 
this  FIR  filter  is  also  compared  with  He(f)  on  Figure  4-2.  The  agreement  is  very 
good  with  an  error  smaller  than  0.2  dB. 


4.2  HILBERT  TRANSFORM 

In  the  previous  section,  a  technique  for  designing  a  filter  matching  a 
given  transfer  function  H(w)  was  presented.  It  relied  on  the  full  definition  of 
H(«) ,  ie.  both  its  real  and  Imaginary  components  (or  alternatively  both  magnitude 
and  phase) .  Filters  are  frequently  defined  in  terms  of  magnitude  or  squared- 
magnitude,  and  in  many  cases,  it  is  difficult  or  not  practical  to  measure  both 
components.  For  instance,  with  the  use  of  a  sine-wave  generator  and  a  power 
meter,  it  is  easy  to  measure  the  magnitude  of  the  response  of  a  digitizer,  but 
not  possible  to  obtain  the  phase  information.  The  phase  response  cannot  be 
simply  ignored  if  a  stable  and  causal  system  is  desired. 

Linear  systems  theory  ([7],  [9],  [15],  [18])  shows  that  for  stable  and 
causal  systems,  the  real  and  imaginary  parts  of  the  Fourier  transform  are  related 
to  each  other  by  the  Hilbert  transform.  In  particular,  for  minimal  phase 
systems^,  the  following  relation  for  H(w)  -  H]|(w)  +  jH3(w)  is  obtained: 


^  A  minimal  phase  system  is  one  for  which  all  the  poles  and  zeros  of  the 
transfer  function  H(z)  are  inside  the  unit  circle. 
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Figure  4-2.  Response  of  filter  for  compensating  cable  losses. 
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Figure  4-3.  Impulse  response  of  filter  to  compensate  cable  losses. 


28 


H,(w)  »  HjCw) 


(29) 


where  H  denotes  the  Hilbert  transform  of  H.  Similarly,  the  relation  between  the 
magnitude  and  phase  can  be  easily  demonstrated  if  we  define  H'(u)  as  the 
logarithm  of  H(w) ,  we  obtain: 

H'(w)  -  ln(  |H(w)|*e^^"’  ) 

'  '  (30) 

=  ln|H(a>)  I  +  j  H^(w) 

which  can  be  related  to  (29). 

As  with  the  Fourier  transform,  a  discrete  Hilbert  transform  (DHT)  can  be 
defined.  However,  it  is  easier  and  more  efficient  to  relate  the  DHT  with  the 
DFT.  The  properties  of  the  Hilbert  transform  are  developed  in  [15]  ,  but  the  most 
important  ones  are  summarized  here: 

The  Hilbert  transform  g(t)  of  a  real  signal  g(t)  is  real. 

The  duality  property  states  that  if  the  Hilbert  transform  of  a  signal  g(t) 
is  g(t),  then  the  Hilbert  transform  of  g(t)  is  -g(t). 

And  more  Important,  the  Fourier  transforms  of  g(t)  and  g(t)  are  related  as: 

G(f)  -  -Jsgn(f)G(f)  (31) 


where  sgn(f)  is  the  slgntun  function.  This  relation  allows  implementation  of  the 
DHT  simply  by  using  the  already  developed  DFT  or  FFT  algorithms. 

For  a  sequence  of  N  points  formed  from  a  signal  and  its  mirrored  conjugate, 
the  slgnum  function  is  defined  as: 


sgn(n) 


1  for  n  -  1,  2,  ....  ^ 

•  -1  for  n  *  .^+1,  ....  N-1 

0  for  n  *  0, 


(32) 
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From  the  properties  presented  above,  a  simple  procedure  can  be  established 
to  recover  the  phase  Information  from  the  magnitude  alone : 

Form  a  N-point  sequence  (x)  of  the  magnitude 
Take  the  natural  logarithm  of  (x) 

Mirror  {x}  with  its  conjugate,  thus  obtain  an  extended  sequence  {x„,} 
of  M-polnts. 

Take  the  FFT  of  {x,j) 

Multiply  with  sgn(n) ,  then  with  -j 

Take  the  real  part  of  the  inverse  FFT 

Truncate  back  to  N-point 

The  example  shown  on  Figure  4-4  illustrates  this  technique.  The  phase  of 
a  resonating  2"'*-order  system  is  reconstructed  from  its  magnitude,  shown  on  top. 
Due  to  the  periodic  nature  of  the  DFT,  the  phase  obtained  converges  back  to  zero 
at  the  maximum  frequency,  but  as  demonstrated  in  the  previous  section,  this  can 
be  interpreted  as  a  fixed  amount  of  delay  which  has  no  effect  on  the  design  of 
a  FIR  filter.  After  compensating  for  that  delay,  the  phase  obtained  with  the 
Hilbert  transform  is  in  good  agreement  with  the  known  phase. 
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Figure  4-4. 


Use  of  the  Hilbert  transform  for  phase  recovery. 
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5.0  COWCLPSIOH 


It  was  shown  that  proper  digital  signal  processing  can  significantly 
improve  the  accuracy  of  the  measurements  or  extend  the  useful  bandwidth  or  range 
of  an  instrument.  Algorithms  were  developed  to  take  into  account  the  transfer 
functions  inherent  to  each  of  the  various  sensors.  Other  techniques  were 
developed  to  compensate  for  errors  present  in  the  measurement  system. 

It  is  necessary  to  integrate  the  output  of  electromagnetic  field  sensors 
responding  to  the  time -derivative  of  the  field.  An  alternative  approach  to 
perform  this  integration  was  developed.  A  specially  made  device  (called  a 
partial  integrator)  was  designed  and  used  in  conjunction  with  digital  filtering. 
Results  showed  that  the  bandwidth  of  the  digitizers  was  extended  by  allowing  the 
use  of  slower  sweep  rates,  and  that  the  dynamic  range  was  extended  by  compressing 
the  peak  of  the  output  and  thus  expanding  the  other  portions  of  the  waveform. 

It  was  shown  that  some  field  sensors  may  be  used  to  measured  the  field 
directly  (and  not  the  time-derivative).  These  measurements  are  affected  by  a  low 
frequency  cut-off  and  an  algorithm  was  developed  to  compensate  this  effect. 

A  technique  for  designing  FIR  filters  with  arbitrary  transfer  function  was 
introduced.  This  technique  is  particularly  well  suited  to  model  the  effect  of 
signal  transmission  through  long  cables  and  also  to  fine  tune  the  response  of  HR 
filters.  In  many  cases,  very  good  response  (error  less  than  a  fraction  of  a  dB) 
can  be  obtained  with  few  filter  coefficients  (between  10  to  50) . 

The  relation  between  the  Hilbert  transform  and  the  real  and  imaginary  parts 
of  a  transfer  function  was  presented.  This  relation  was  extended  to  relate  both 
the  magnitude  and  phase  of  the  function,  thus  allowing  to  reconstruct  the  phase 
information  from  the  magnitude  alone.  This  relation  allows  the  design  of  a  FIR 
filter  to  model  a  system  for  which  only  the  magnitude  can  be  measured.  A 
practical  implementation  of  the  discrete  Hilbert  transform  was  also  presented. 
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